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ABSTRACT 

We obtained spectra, in the wavelength range A = 995 — 1769 nm, of all four known planets orbiting 
the star HR 8799. Using the suite of instrumentation known as Project 1640 on the Palomar 5-m 
Hale Telescope, we acquired data at two epochs. This allowed for multiple imaging detections of the 
companions and multiple extractions of low-resolution {R ^ 35) spectra. Data reduction employed two 
different methods of speckle suppression and spectrum extraction, both yielding results that agree. 
The spectra do not directly correspond to those of any known objects, although similarities with 
L and T-dwarfs are present, as well as some characteristics similar to planets such as Saturn. We 
tentatively identify the presence of CH4 along with NH3 and/or C2H2, and possibly CO2 or HCN in 
varying amounts in each component of the system. Other studies suggested red colors for these faint 
companions, and our data confirm those observations. Cloudy models, based on previous photometric 
observations, may provide the best explanation for the new data presented here. Notable in our data 
is that these presumably co-eval objects of similar luminosity have significantly different spectra; the 
diversity of planets may be greater than previously thought. The techniques and methods employed 
in this paper represent a new capability to observe and rapidly characterize exoplanetary systems in a 
routine manner over a broad range of planet masses and separations. These are the first simultaneous 
spectroscopic observations of multiple planets in a planetary system other than our own. 
Subject headings: planetary systems; planets: individual (HR 8799bcde); stars: individual (HR 8799); 

instrumentation: adaptive optics, spectrographs, coronagraphs; methods: data 

analysis; techniques: spectroscopic 



1. motivation and properties of HR 8799 

The star HR 8799 (HD 218396; V432 Pegasi), an A5V 
7 Dor- type variable star ( Rodriguez fc Zerbi||1995 Aerts 
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et al.||1998| Gray fc Kaye||1999D at 39.4pc ( |van Leeuwen 
2UU7|), has been the target of inany observations in the 



past few years due to the discovery of several fa i nt com - 
panions in orbit about it (Marois et al. 2008 2010b). 
The nature of these sources — whether they are planets 
or brown dwarfs (e.g Moro-Martm et al.|2010 ) — has been 
subject to significant debate. Various observations have 
attempted to constrain the age of the system (and thus 
the masses of the companions based on theoretical cool- 
ing models); to obtain as many uncontaminated photo- 
metric measurements of the companions as possible; and 
to determine their orbits. Here we refer to them as plan- 
ets while acknowledging that some researchers still have 
doubts about the use of this label. 

In order to understand the physics and chemistry 
of these objects, spectroscopy is necessary. Herein we 
present the results of high- contrast, direct, spectroscopic 
imaging observations of this exosolar system after a short 
review of its properties, as observed to-d ate. A compan- 
ion paper (Paper H; Pueyo et al. 2013) reports the as- 
trometric observations of the system, and a third paper 
(Veicht et al. 2013) examines the dynamical stability of 
the system. 

For reference, detailed reviews of hig h- contrast, coron- 
agraphic observations can be found in |Traub fc OppeiT 
heimer (2011), Absil & Mawet (2010] ) and Oppenheimer 

^rnhkieyii^ooy " 



1.1. Variability and Rotational Velocity 



Oppenheimer et al. 



TABLE 1 
Observations 



Date 
(UT) 


Julian 
Date (start) 


Hexp X texp 

(s) 


Estimated 
seeing (")^ 


Estimated 
Strehl^ 


CAL 
Mode^ 


rms-WFE 
(nm)d 


Planets 
Detected 


14 June 2012 


2456092.972859 


5 X 549.91 


1.4 


65% 


<l> only 


5 


bcde 


15 June 2012 


2456093.971030 


5 X 366.61 


1.1 


67% 


<l> only 


4 


bcde 


5 October 2012 


2456205.706470 


18 X 549.91 


1.8 


44% 


<l> and ^ 


14 


c 



^ Seeing was estimated using short images with no AO correction and is quoted as the FWHM of the PSF as measured 
at A = 1653 ± 12.5 nm. 

^ Strehl was estimated at A = 1653 it 12.5 nm using unocculted, AO-corrected images of the primary star taken just 
before coronagraphic occultation and long integrations and does not include the effects of CAL on the mid-spatial 
frequencies in the image. 

c "<|) only" means phase correction only, while means full electric field conjugation was applied for optical speckle 

control. 

^ rms wave front error in spatial frequencies from 5 to 32X/D as measured by CAL, without seeing. 



The class ification of HR 8799 as a 7 Dor variable (Kaye 
itl[L|[l999 ') was secure by 1997. Periodicity of the varia- 



t ions from Hipparcos data ( Aerts et al.][l998 ) and a long 
photometric campaign revealed multiple periodic signals 
with frequencies ranging from 0.2479 to 1.9791 days with 
amplitud es up to 1-5% an d a range of ±4% in the longest 



mode (Zerbi et al. 1999). The derived stellar proper- 
ties from thi s study are L = 5. QI/(t^, Teff = 7230K and 
R = I.URq ( jBaines et al1[2QT2l ). 

Large time-series data sets of photometry and spectra 
have been taken on HR 8799 and a number of asteroseis- 
mology analyses have been conducted, first to constrain 



the age, unf ortunately with only weak constraints, (Moya 



et al. 



2010) and later to constrain t he inclination an^le 



of the star s rotation axi s {i > 40°) ( [Reidemeister et aL 



2009[ [Wright et al.||201lD . The rotational velocity ol the 



star is vsini = 37.5 to 49 km s~^, based on a variety 
of different m easurements ( Zorec & Royer||2012 Kaye & 
' Strassmeier|p'998 J . 



1.2. Age and Planet Masses 

Though still under debate, several indic ators suggest 
that the star is y ounger than 100 Myr (Marois et al. 
2010bl iDoyon et al. 2010, Zuckerman et al. 2011. Malo 



et al.||2013||Moya et al.||2Q10p . One estimate" of t his sys 
t em's age comes from two indepen dent studies ( [Doyon 



Zuckerman et al.||2011| ), each of which used 



galactic space motions and spectroscopic age indicators 
to conclude that HR 8799 has a very high probability 
(~98%) of being a me mber of the 30 Myr Columba As- 
sociation. In contrast, Hinz et al. (2010) convincingly 



demonstrate that this association is poorly supported 
and that the age of HR 8799 cannot reliably be asso- 
ciated with Columba. Earlier works argued for a system 



age youn ger than 50 -60 Myr (Mar ois et al.l 
[meister et al. 2009). While si gnificantly d ; 
the 1 Gyr estimate put forth by Moya et al 



2008 Reide- 



itterent from 



Gyr estimate put for th bylMoya et al. | (|2010[) using 
asteroseismology data, the Moya et al. ( [2010 ) work does 
not explicitly exclude such a young age. 

The age of the system is important because it con- 
strains the masses of the planets when used in con- 
junction with theoretical cooling models and the mea- 
sured luminosit y of the objects . It is interesting to note, 
as discussed in Skemer et al. (2012), that some of the 
strongest constraints on the age 01 the HR 8799 sys- 
tem come from the dynamical simulations of the four 



planets, combined with their broadband photometry. As 
Gozdziewski & Migaszewski (2009) discuss, only a small 
fraction of three- body systems resembling HR 8799 can 
remain stable longer than ~100 Myr. Further, the stabil- 
ity of a four planet system for ~10-100 Myr is consistent 
with lower masses for the companion s ( Gozdzi ewski 
Migaszewski||20091 [Reidemeister et aL][2Q09, .Fabrycky 



Murray-Clayj|2010( [Moro-Martm et al [j201Q[). Specifi- 
cally, t he study presented in 'G ozdziewski ^ Migasze wski 
(|2009[) finds only very limited dynamical configurations 
tnat foster long term stability: such stability will occur 
if the planets exhibit low-order, two- or three-body mean 
motion resonances. 

These conclusions are consistent with the work pre- 
sented by Fabrycky & Murray-Clay (2010), which find a 
longer system lifetime if the two innermost planets main- 
tain a 2:1 commensurability. Further, the persistence of 
this 2:1 resonance to the current epoch would put a spe- 
cific constraint on the masses to be less than 10 Mj^p. 
Taking this one step further, [Fabrycky fc Murray- C lay [ 
(2010) point out that these dynamical constraints serve 
as early te sts of the so called "hot start" models (Baraffej 
et al.[[2003). Specfically, if all three planets exhibit low- 



order commensurabilities with each other, the hot start 
models are easily consistent with the measured luminosi- 
ties and inferred masses. If, on the other hand, only 
the two innermost planets have a 2:1 commensurability, 
the hot start models are only barely applicable. These 
studies point to the promise of using high contrast imag- 
ing to place empirical constraints on evoluti onary models 
for th e brightness of substellar objects (e.g. Crepp et al. 
2012[) . Nonetheless, the p ublished absolute i i-band mag- 
nitudes, Mh ^ 14 to 15 ( Marois et al.[[2008 ), are consis- 
tent with standard coolin g evolutionary models for 30 
Myr ( Burrows et~aL][1997[ ). 
Further, the d egree of excess infrared emission mea- 



sured at 24 jam (Su et al. 2009) is consistent with an 
age significantly younger than 1 Gyr. Observations of A- 
stars with circumstellar material plus well-constrained 
ages show a marked decline with age in the strength of 
the 24 jam infrared excess (e.g. |Rieke et al.n2005 ). This 
empirical res ult can be well- mode lled by dynamical sim- 
ulations (e.g. Wyatt et al. 2007), illustrating the clear 
decay of this hot dust emission due to the shorter dynam- 
ical time scales associated with the innermost portions of 
these systems. 
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Fig. 1. — Images demonstrating the optical suppression of speckles through wave front phase error sensing and control. Left: coronagraphic 
image of HR 8799 from 15 June 2012 without the CAL system in operation. Right: the CAL system iterated 5 times with phase conjugation 
only ("<l>-only") to suppress speckles within the 32X/D distance from the star (within the blue square). In the image to the right, the 
control region of PALM-3000 (32A/D) is clearly apparent as the darker square area centered on the star and extending to the edge of 
the control region (the blue square), where the speckles brighten again. Images are both 549-s long and are single slices at A = 1224 nm. 
Color scale is logarithmic and identical in both images. North is up and East to the left, with an image width of 3(^8. The planets, whose 
locations are shown in Fig.|2] have brightnesses of 3.2, 3.3, 2.9 and 3.7% of the mean speckle brightness in the vicinity of each planet, while 
the CAL system is operating, for b through e respectively. Thus, advanced image processing techniques are required to study them. 



Based on the ages of 30, 60 and 430 Myr the four com- 
panions have masses in the fohowing ranges: 7, 13 and 
34 Mjup for c, d and e; a nd 5, 14 and 23 M jun for b 
using the coohng models of Baraffe et al. (2003). An im- 
portant cave at exists regarding the age of the system, as 



ig tne i 

discussed in |Moya et al.| ([20107. If the inclination angle 
of the star is near 50"" and the mass of the star is 1.45 
M0, the age range from asteroseismology analysis is 26 
to 430 Myr, with the larger ages pushing the masses of 
all components into the "b rown dwarf" mass range. As 



|Moro-Martm et ahj (2010) point out, the inclination of 
the star is probal5Iy smaller than the value initially as- 
sumed. Nonetheless, spectroscopy can also be used to 
constrain mass, or more accurately, log (g), as described 



et al. (2011). 



in, for exan iple, Marley et al. (2012) and [Madhusudhan 



1.3. Previous Spectroscopic Measurements 

A number of studies have accomplished spectroscopy 
of the b component in the L, K and H ba nds (Janson 
et al.l|2010[ [Bowler et aT]|2010[ [Barman et al.[|20TT| nn 
addition, preliminary observa tions of the ccomponen t 



in K band were present ed by jKonopacky et al 
In 



Janson et al. 



(2010 



( |2013[ ). 

a noisy L-band spectrum of 
ultimately sub-sampled to very low resolution, sug- 
gested that either improvements in the treatment of dust 
are required to model the object properly, or significant 
non-equilibrium ch emistry is present in b's atmosphere. 
Bowler et al. (2010j) ob tained a K-band spec trum of b, 
and, combined with the Barman et al. (2011) da ta in H 
and K, a n umber of modeling papers, notably "Marley' 
et al. (2012), suggest that complex cloud structure, high 



metallicity and/or non-equilibrium chemistry are needed 
to understand this object. In any case, b shows little 
or no methane absorption, possibly some CO and the 
shape of the ij'- band peak suggested low-gravity ( [Bar-] 
man et al. 2011). Unfortun ately at the time of wri tmg, 
the c component spectrum (Konopacky et al. 2013) and 
its analysis were not available in any publication. 

These previously published measurements are all of the 
outer most planet in the system. However, this currently 
unique system has four directly detected planets. Thus it 
is ideal for the study of multiple planets around a single 
single star to investigate the diversity of planets. This is 
the purpose of this paper. 

1.4. Orbital Characteristics 

HR 8799 has become a benchmark target of observa- 
tion in the field of direct exoplanet imaging to demon- 
strate whether a given project or technique is capable 
of finding even fainter objects, relative to the star, than 
the four planets orbiting this star. As a result numer- 
ous images of the system have been obtained main! 
the near IR through L- and M-ba nds (Hinz et al. 2010 



Esposito et al. 2013 



2011a 
et alJ 



Currie et al 



Currie et aL]|2011| [Hin kley et 



20121 IG alicher et al.||2Ql In Janson 



20101 ILafreni^re et al. 2009; Marois et al.||2008 



Serabyn et al. '2010 ; Ske mer et al. [20121 [Soummer 



2010b . 

et al. 2011; Currie et al. 2012). The conclusion of this 
work from various different instruments and telescopes 
indicates that the three outer planets are roughly copla- 
nar, orbiting at an inclination angle of i ^ 27.3 to 31.4°, 
with eccentricities below 0.1 and semi-major axe s of 68, 
42 and 27 AU for b, c, and d respectively (see esp. Soum- 
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mer et al. 2011 Currie et al. 2012). Orbital motion of 
e is poorly constrained at present, and only two groups 
attempted dynamica l studies including a ll four planets 
dEsposito et aL]|2013( [Marois et aLl|2010bD . 

I'he astrometry oi this system is complicated by the 
background speckle field and, due to the extensive na- 
ture of the data reduction , we present that analysis in 
a second paper (Paper II; Pueyo et al. 2013 in prep.). 
However, HR 8799 will likely be the first exoplanetary 
system with a complete characterization including dy- 
namical masses from astrometry through atmospheric 
analysis via spectroscopy, part of the purpose of com- 
bining imaging with spectroscopy in a single instrument 
such as Project 1640. We defer further discussion of the 
orbital characteristics to Paper II. 

1.5. Circumstellar Material 



The HR 87 99 system was first identified by "Zucker-' 
man fc Song| (2004 ) as a source with a prominent IRAS 
infrared excess, indicating the presence of cool dust. 
Shortly thereafter , the analysis of Sp itzer spectroscopic 
data presented in Chen et al. (2006) indicated that the 
debris structure in this system was organized into two 
spatially distinct components. After the initial discovery 
of the planets in the system, Su et al. (2009) presented 
the analysis of much deeper Spitzer spectroscopy as well 
as photometry at 24 and 70 jam. In addition to con- 
straining the spatial extent of the dust belts, this work 
identified an extensive outer halo of dust particles at hun- 
dreds of Alls, and these larger regions were later probed 
thr ough observations pres ented in Hughes et al. (2011), 
and Patience et al. (2011) at 350 and 880 /xm, respectiv- 
ley. 'I'hese observations were used to resolve the outer 
debris structure surrou nding the planets, solidifying the 
framework presented in'Su et al. ( |2009| ). A coherent pic- 
ture of the overall architecture of the system near the 
time of the planet discovery, including the third outer 



component to the d ebris disk, was synthesized in |Reide- 
meister et al.|(|2009|) 



I'he inner structure of the HR 8799 debris disk is 
shaped through dynamical interactions with the four 
planets. Specifically, the radial gap in the debris struc- 
ture from 15-90 AU is clearly consistent with clearing 
caused by the four known planets. However, the source 
of clearing interior to 6 A U in the inner warm b elt re- 
mains a mystery, with only Hinkley et al.| (2011a) plac- 
ing firm constraints in this region. It a p laneta ry mass 
companion is responsible for the clearing, [Hinkl ey et al. 
( |2011a[ ) put a firm upper limit of ~11 Mjup on its mass, 
rulmg out the possibility of a companion star or brown 
dwarf between 0.8 and 10 AU. Spectroscopy in the pub- 
lished literature has not revealed any evidence of binarity 
either. 

2. INSTRUMENTATION: PROJECT 1640 



Project 1640 is desc r ibed in Oppenheimer et al. ( 2012|) , 



Hinkley et al. 
at the level o: 



( |2011cD , [Hinkley eta l. 12008) and in detail 
circuit diagrams, cryogenics, control soft- 



ware, interfaces and opto- mechanical design in |Hinkley| 
( 2Q09|). The latest system per formance metrics are given 
in [oppenheimer et al. (2012), including on-sky contrast 
measurements. These are described in relation to other 
projects in high-contrast imaging in Mawet et al.,(2012^ . 



in particular, their Fig. 1. In summary, the system is ca- 
pable of producing images with a speckle fioor at roughly 
10~^ at separation from a bright star (or 10~^ in 
the lab). This is achieved through the coordinated op- 
eration of four optical instruments: a dual deformable 
mirror, adaptive optics (AO) system with 3629 a ctively- 



controlled actuators, called PALM-3000 (Dekany et al. 



2007 2006, and 2013, in pr eparation); an apodized p upil. 



Lyot coronagraph (APLC; Soummer et a l. 2009; Si vara- 
"makrishnan fc Lloyd||2005[ |!Soummer| |2005( (Soummer 
et al. 2003a b), the design details of which are given 
in ( [Hin kley 2009); a Mach-Zehnder interferometer that 
senses and calibrates, through feedback to PALM-3000, 
residual path-length and amplitude errors in the stel- 
lar wave front at the coronagraphic occulting spot for 
optimal diffr active rejection of the pr imary star's light 
(GAL; Vasisht et al. 2012; Zhai et al. 2012); and an in- 
tegral field spectrograph that takes 32 simultaneous im- 
ages with a field of view of 3^8 x 3^8 spanning the range 
A = 995 - 1769 nm with a bandwidt h of AA 24 9 nm 
per imag e (IPS, Fig. [l| [Hinkley et a L][2011c[ [20081 [Hink- 



ley||2009| [Oppenheimer et al.||2012f 7 Aside from tech 



meal advances m high-contrast imaging, numerous re- 
sults from the project include, among others, the discov- 
ery and astrometric and spectroscopic characteriz ation 



of the Alcor AB sys tem (IZiinmerman et al. 2010), the 
a Ophiucus system (Hinkley e t al.[25llb[ ), tne (" Virgi- 
nis companion (Hinkley et al. 2010[), and compre hensive 



spectral studies of the compani on of FU Orion is ( Pueyo 
et al.^ 2012b) and Z CMa (Hink ley et al.|[2QT3 |. 

Raw science data generated by Project 1640 are in 
the form of 2040 x 2040 pixel images containing 37146 
closely packed spectra roughly 30.4 x 3.2 pixels in extent. 
These images are processed into data cubes with dimen- 
sions R»A. 
( [2QlT| 



5 and A, as described in [Zimmerman et al. 
Cross talk between adjacent spatial pixels and 
across wavelengths is less than 0.4% flux contamination 
on average with a maximum of 1.4% in the water band 
between J and i^-band, as confirmed using monochro- 
matic sources and testing of the cube extraction soft- 
ware (see Oppenheimer et al. 2012). This level of cross 
talk does not introduce systematic errors in spectral ex- 
tractions or image properties that would affect relative 
photometry at a level higher than 2% to be conservative. 

Because the planets we observed have relative bright- 
nesses of 3.2, 3.3, 2.9 and 3.7% of the mean speckle 
brightness in the vicinity of each planet, for b through e 
respectively, and the speckles vary spatially and tempo- 
rally, advanced image processing algorithms are required 
in addition to the complex suite of hardware and software 
described above (Figs. [2] and [3|. 

3. OBSERVATIONS AND DATA REDUCTION 

Target acquisition for Project 1640 involves placing the 
star within the field of of view with telescope pointing, 
optimizing the AO correction, obtaining images of the 
star when it is not occulted by the coronagraphic spot 
(called Core images — acquired only if the star does not 
saturate the detector in short exposures), acquiring the 
star behind the coronagraphic spot (using a fine guidance 
sensor that centers the star on the spot), enga ging the 



CAL system to dim quasi-static speckles (e.g. Hinkley 
et al.|[2007) and then acquiring deep exposures sensitive 



enough to detect the speckle fioor. Core images are used 
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KLIP Algorithm 





. > , r ."^ ■■ - ^ J. ' 



Fig. 2. — Speckle suppressed images of the HR 8799 system from 14, 15 June 2012. Left: Speckle suppression by the KLIP algorithm 
for June 2012 (summed over both dates) in /f-band. Right: same as left but speckle suppression was achieved by the S4 algorithm. These 
images both represent a total of 4582. 61-s of exposure time. In October 2012 only the c component was detected due to poor seeing 
conditions and an incorrect field-of-view setting that vignetted the position of the b component (see text). Because of differences in the 
two algorithms, the blocked out areas in the center of the images are different. 

for relative photometry and spectral calibration. If Core 
images saturate, a fainter star from the IRTF spectral 
library (Rayner et al. 2009) is used for these calibrations. 
In the case ot HK 8799, the Core image does not saturate 
in short (1.45-s) exposures. 

The CAL system senses the full electric field of the inci- 
dent light and can operate in two modes by feeding back 
different control signals to PALM-3000: (1) suppressing 
speckles only due to wave front phase errors in Table 
IT]), which generates a dark square region centered on the 
primary star with an angular size of 32X/D exten ding 
from the star in all fo ur cardinal directions (Fig. [ij |0p- 
penheimer et al.|2012 ), where D is the telescope diameter 



as modified by any pupil-plane stops in the optical sys- 
tem; or (2) suppressing speckles due to both phase and 
amplitude errors in the wave front to produce a much 
darker region in only half of the field of view, extending 
from the star to 32 X/D but only in three of the ca rdinal 
directions (NSEW; Vasisht et al.||2012| |Zhai et al.|[2012( ). 
While the second mode works, the observations used in 
this paper were primarily with phase-only ($) control, 
with about half of the October 2012 observations using 
electric field conjugation. 

Between long exposures, short exposures (tint = 30- 
s) for the purposes of p recision astrometry {a o:^ 2 mas; 
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Zimmerman et al. 2010p were obtained in which a 150 



nm amplitude sine wave was induced on the high-order 
deformable mirror in two orthogonal directions with 22 
cycles each over the telescope pupil. These create four 
spots, at a separation of 22\/D^ that permit precise de- 
termination of the position of the star behind the coro- 
nagraphic mask. The details and use of these spots 
have been described else where and will be dealt with in 
more detail in Paper II (^Digby et al.||2006l Marois et al. 



Fig. 3. — S4 detection map for June 2012. Orientation is N up 
and E left. Coordinates are in pixels and the detection map is a vi- 
sualization of the correlation between a companion PSF, obtained 
experimentally, and the residual image (Fig. 12] right). Compan- 
ions are found with higher fidelity using the oetection map than 
the residual image alone, because it rejects many apparent false 
positives in the residual image. 



2006[ Sivaramakrishnan ^ Oppenheime'r immer- 



man et al.|[2QTQ ). 

Observations of HR 8799 on 14 and 15 June 2012 com- 
prised a total of 46 and 31 minutes of exposure time, 
while 165 minutes of integration time were obtained on 
5 October 2012, respectively. Details of the observations 
and conditions are provided in Table [l] Astrometry im- 
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ages (as described above) were not obtained in June 2012 
but were taken in October 2012. Examples of one wave- 
length slice of two different cubes are shown in Fig. [T] be- 
fore any speckle suppression post-processing. The right 
image in Fig. [l] shows the effect of the GAL system, in 
<l>-only mode, on the stellar point spread function (PSF), 
creating a darker square centered on the star, within the 
AO control region (32 A/D). 

Our experience with the PALM-3000 AO system shows 
that Project 1640 can only achieve its full sensitivity 
when the natural seeing is below about 1^7, which was 
the case during the June 2012 observations. In Octo- 
ber 2012 however, conditions were considerably worse. 
In addition, due to a slight misalignment between the 
IFU and coronagraph, the b component was outside the 
field of view in October 2012. As a result, the October 
2012 data reveal only the c component of the HR 8799 
system. It is important to note that the inner planets 
were not seen in October, despite the significantly longer 
exposure time, because of the larger residual wave front 
error (WFE) of 14 nm rms as opposed to the 4 and 5 
nm achieved in June. Speckle brightness is a very strong 
function of rms- WFE when it is well below the 50 to 150 
nm achieved by most AO systems at these wavelengths. 
As such, the speckles were more than one hundred times 
brighter than the planets in October, making them invis- 
ible, even to the advanced speckle suppression techniques 
described below. 

4. SPECKLE SUPPRESSION 

Raw images (Fig. [T]) from the extracted cubes do not 
directly reveal the companions of HR 8799. The speck- 
les in the outer parts of the PSF have to be removed in 
an effective manner or modeled precisely. We have used 
two different techniques, each of which employs a new ap- 
plication of principal component anal ysis in astronomy . 
The first method, called KLIP (S oummer et al.||2012^ , 
operates to suppress speckles within a given cube, re- 
moving the speckles from each slice individually. Details 
on this algorithm are fflven i n Ap pendix [A] The second 
technique, S4 (Fergus et al. 2013), uses the full depth 
of the diversity of measurements in the data, including 
modeling the changes in the speckles in time, from cube 
to cube. See Appendix [B] for more details. 

Fig. [2] presents images of the system after speckle sup- 
pression with KLIP (left) and S4 (right). Note that all 
four companions are visible, though the S4 algorithm ac- 
complishes higher signal-to-noise ratio (SNR). For each 
component the maximum SNR at any wavelength is 15.4, 
29.1, 9.7 and 8.9-cr for b, c, d and e respectively (Table p|. 

The S4 algorithm uses the residual image to identify 
possible companions in a given set of data automatically. 
This is done by cross-correlating the residual images with 
a model companion PSF, which is actually real data ac- 
quired with the laboratory white-light source over a grid 
of 25 positions on the field of view. These 25 PSFs are 
used to account for any variation in PSF shape and size 
over the range of positions and wavelengths in the cubes. 
The algorithm assumes a "white" spectrum and gener- 
ates, via the PSF cross-correlation, what we call the "de- 
tection map." This map for the June 2012 observations is 
shown in Fig. ^ (for 200 components in the PGA model; 
see Appendix Pi . In this image, positive correlations can 
be as sharp as a single pixel, and one can see that all 



four planets are clearly detected, with higher SNR than 
in the residual image. This is partly why S4 is such an 
effective speckle suppressor. 

Although no a priori information on the locations of 
the four components was used in the detection algo- 
rithms, the e component would not have been discovered 
with these techniques, at least on simple examination of 
the detection map and the residual image produced by 
S4. However, there is a priori information on this system 
from other observations, so we were able to retrieve the e 
component. We note that an extension to the algorithm 
examines the extracted spectra of all peaks in the de- 
tection map, comparing them to the general background 
spectrum to validate whether a detection is real. If this 
were implemented at this point, e would certainly have 
been flagged and discovered without data from other ob- 
servations, because its spectrum is markedly different 
from the background spectra at the same angular radius 
from the star (see Another possible way to increase 
detection probability would be to use an assumed input 
spectrum, other than white, for the type of companion 
being sought. These improvements are discussed in |Fer- 
,gus et al.. ( ,2013J . 

5. PHOTOMETRY 

From the cubes of speckle suppressed data, photometry 
can be extracted. KLIP data have been shown to co ntain 
negligible bias in photometry (jSoummer et al.|2012|) , un- 
like the LOCI algorithm which was used for previous de- 
terminations of photome try of this system (Maro is et al.| 
2000| |2003| |2008[ |2010aD . Thus it is the best choice tor 
deriving photometry here. The data are in units of con- 
trast, so they can easily be converted into photometric 
measurements using th e absolute m agnitude of the pri- 
mary star (Mh = 2.30; Cutri et al.| ^2003). We estimate 
the broadband photometry by propagating these spec- 
tra through the J and H MKO filter bandpasses. For 
the channels where a detection cannot be established we 
the value of the spectrum is not included in the pho- 
tometry. Our findings are summarized in Table [2] and 
are in good agreement with published values, albeit with 
slightl y lower J-band flux than reported in Marois et al.| 
(2008) for HR 8799 c. Interestingly, this makes c closer 
in color to b and d. The d and e components were not 
detected in the J-band. 

We note that all data from each epoch were needed 
to make these measurements and there is no statistically 
significant difference between the photometry for c in 
June and October. 

6. SPECTROSCOPY 

We derived spectra of each of the four point sources 
using both KLIP and S4, with neither algorithm making 
any assumptions of spectral shape. Both algorithms yield 
consistent results for the brightest detections in the wave- 
length range. These are the critical wavelength channels 
for the photometry presented above, meaning that either 
algorithm could be used for photometry, but due to the 
calibration of the KLIP output, it was more straightfor- 
ward to use KLIP data. For the fainter channels, the S4 
algorithm is superior at obtaining the spectra, because it 
uses all of the information in the original data cubes — 
KLIP treats the channels independently and thus does a 
poorer job of speckle suppression. For these reasons we 
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TABLE 2 

Absolute Photometry of the HR 8799 System 





HR 8799 b 


c 


d 


e 


Project 1640 Mj 
Published^ Mj 


16.48 ± 0.18 
16.30 ± 0.16 


15.36 ± 0.21 
14.65 ± 0.17 


> 15.5 
15.26 ± 0.44 


> 13.2 
No reported value 


Project 1640 Mr 
Published" Mr 


15.11 ± 0.12 
14.87 ± 0.17 


14.20 ± 0.15 
13.93 ± 0.17 


13.69 ± 0.19 
13.86 ± 0.22 


13.30 ± 0.27 
13.53 ± 0.44 



^ Published val ues are from|Marois et al.| ( |2QQ8| ) and |Marois et ar| ( |2QlQb| and the primary 
has Mh = 2.30 ( jCutri et al.|2003| ). 



decided to conduct spectral analysis on the S4 extrac- 
tions, which are presented in Fig. |4]in arbitrary, normal- 
ized units of Ta, which is commonly done for comparative 
spectroscopy. The spectrum of the c component from the 
October 2012 observations is also shown in Fig. [4j It is 
consistent with the June 2012 spectrum in the JY-band 
within the error bars, except for a minor 2-cr difference 
in a new feature in the points flanking 1600 nm (see 
We are unaware of observations that would correlate tms 
difference with variability in the primary star, such as a 
large increase in UV flux that might induce photochem- 
istry to produce a new absorption feature. In the October 
2012 data the c component was only marginally detected 
in the J-band due to the poor seeing conditions. Thus 
the apparent, but not real, discrepancy between June 
2012 and October 2012 at roughly 1300 nm. Ah other 
points are consistent within the errors. For components d 
and e, we have excluded the four points obtained within 
the telluric water band between 1350 and 1430 nm, be- 
cause what little astronomical signal is present is dom- 
inated by starlight. For b and c this was not an issue 
because of their much more distant separation from the 
star. The other channels had negligible cross-talk from 
the star, although the points short ward of 1350 nm for 
e and d are considered upper limits in our analysis. 

As indicated in Appendix [Bj the exact range of prin- 
cipal components and size oiregions used depends upon 
seeing and other weather conditions for optimal SNR. 
The extractions for October required a larger number of 
principal components. (See Appendix [B|) 

7. TESTS OF SPECTRAL EXTRACTION FIDELITY 
Because these objects have only be en studied spectro- 



scopically in a limited m anner (see [1.3 [Bowler et al. 
2010 Barman et al.||2011 ) there is only one part of the b 
component's spectrum to compare with to demonstrate 
consistency. In addition, the techniques we used to de- 
tect and measure the spectra are new. Therefore, we 
conducted several tests to provide additional verifica- 
tion that the spectra are not contaminated by some as- 
pect of the instrumentation or the algorithms employed. 
We note, however, that Project 1640 has published nu- 
merous spectra of much brighter companions (ones that 
do not require speckle suppression to detect) and these 
are all consistent with othe r methods of near-IR spec - 
troscopy. See, for exa m] immerman et al.l (12010') 



Hinkley et a l. 
let al.l (|2012[), 



nple 

( 2011b|), [m nkley et al. (2010), Roberts 
al.l 



trum of Titan agrees with other published spectra in 
this wavelength range, exhibiting strong methane and 
water absorption. The fact that bright companion spec- 
tra are reliably extracted does not mean that the same 
techniques can be applied to faint companions below the 
speckle floor. Thus we present four separate tests below 
of the S4 spectral extraction method (Appendix [b|), and 
some statistical comparisons of the four spectra. 

We also note that because both KLIP and S4 repro- 
duce the same features within the errors (Fig. Al ) in the 
spectra and are different techniques, the confidence in 
the spectral extractions is already high. 

The statistical significance of the detections, in the 
form of an average over all A for b and c and just H-hand 
for d and e, and a maximum, both in units of signal-to- 
noise ratio, using the background spectrum as the base- 
line noise against the error bars, is listed in Table [3j We 
note, however, that our methods of determining the er- 
ror bars and using the background spectrum as the noise 
level are conservative, intentionally, leading to over esti- 
mates of the err ors and und erestimates of the SNR (see 
Appendices and Fergus et al. 2013). This statement is 
bolstered by the data shown in Figs. [Sj [8|and|Al 



7.1. Comparison with Other Studies 

The only companion object to HR 8799 with published 
spectral informat ion in this study 's wavelength range is 
component b (Bar man et al.|2011[ ). Our s pectrum in the 
same wavelength region agrees with the [Barman et al. 



(2011) spectrum, except for two points which deviate 
by roughly 2-(j. Fig. [5] shows both spectra overlaid to 
demonstrate the consistency of the two results. The 
current spectrum reveals two distinct features however. 
These are discussed in detail in ^8.3 



Fueyo et all (|2012b) and Hinkley et al. 



( [2013 ). In addition a Project 1640 disc-integrated spec- 



7.2. Other Locations in the Detection Map 

To evaluate the S4 spectral extraction fidelity, for each 
of the four components, we also extracted spectra at 6 
different, randomly selected locations at the same angu- 
lar radius from the central star. These spectra are essen- 
tially fiat and do not reveal the detection of any object. 
Rather they indicate the background against which the 
planets themselves are being detected. These are shown 
in Fig. [4] by the dashed lines for each component. Note 
that in the case of d and e, these background spectra 
show that the components are only barely detected in 
the Y- and J-bands. However, integration of these chan- 
nels permits the photometry for b and c and the upper 
limits for d and e. 
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Fig. 4. — Spectra of the b, c, d and e components using the S4 algorithm. 1-a error bars are indicated on either side of each point. Spectra 
are shown for comparative purposes in normahzed f^. The dotted spectrum of component c is extracted from the October 2012 data and 
shows poorer detection in the J-band due to seeing conditions. All other spectra are from the June 2012 epoch. Black dashed lines are 
spectra of the background extracted at the same radial distance from A for each component. They consist of the average of 6 different 
randomly-selected locations in the azimuthal direction; they represent comparative "noise level" estimates. Tentative identifi catio n of some 
molecular features (excluding water) are indicated at the top of the plot. CO2 is listed with a question mark as explained in ^8.3| and could 
also be attributable to HCN. 



Along these lines, we also derived spectra of fourteen 
other similarly bright peaks in the detection map shown 
in Fig. [3j to see whether they were real sources. All of 
them, which are either not point sources or of lower sig- 
nificance than the bona-fide companions, either showed 
peaks in the insensitive water-band between 1350 and 
1430 nm — indicating spurious starlight contamination — 
or fiat, featureless spectra at levels similar to the nearby 
background, indicating that they are not real sources. 



7.3. Fake Source Spectrum Retrieval 

As another test of the ability of our data analysis to re- 
trieve spectra with high fidelity we inserted multiple fake 
companion sources into the same data cubes we analyzed 
for this paper. The sources were placed at the same over- 
all brightness as each component and at the same radii 
from the star. The fake sources were given the T4.5 
spectrum of the standard 2MA SS J0559-14Q4, chosen for 
its strong molecular features fBurgasser & Kirkpatrick" 
[2006^ , which are useful to understand whether the fea- 



tures are correctly reproduced. The spectrum was resam- 
pled at Project 1640's resolution and wavelength range 
and applied to a cube of a fiducial PSF derived from an 
unocculted observation of the internal calibration white- 
light source. The instrumental response of the system 
was also applied to this PSF and then it was reduced 
in intensity by a factor sufficient to make it of the same 
intensity as each of the b, c, d and e components. We 
then ran this composite data cube through the S4 algo- 
rithm, detected each object with the same significance as 
the real components and extracted their spectra as de- 
scribed in Appendix [B) This was repeated 5 times with 
different, randomly selected azimuthal locations for the 
fake objects at the same radii as the planets. 

The spectra of fake sources extracted matched the in- 
put spectrum, with an rms error of less than 1-a for 
all four cases. All four spectra of the fake sources are 
shown in Fig. |6] In general the extractions are identical 
in shape to the input spectrum and in all cases repro- 
duce small kinks in the input T4.5 object. For reference. 
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TABLE 3 

Average and Maximum Signal-to-Noise 
Ratios^ 





b 


c 

June 


c 

October 


d 


e 


Average 
Max. 


7.0 
15.4 


7.9 
29.1 


4.2 
9.5 


5.1 
9.7 


4.3 
8.9 



^ SNR values use the background and the error 
bars. For the last three columns only points in 
the if-band are included. 




Fig. 5. — Spectrum of HR 8799b. Solid line is our work and 
the dashed line is from Barman et al. (2011). The two spectra are 
consistent within the error bars, except for two points at A = 1538 
and 1567 nm, which deviate by 2.1 and 2.8-cr. The spectra are 
both shown in normalized units of ix. 



each fake component was set to the same 2 to 3% of 
the local speckle brightness as the real planets are (^2j 
i.e. the companions are roughy 40 times fainter than trie 
local speckle background). The average deviation from 
the input spectrum, over all wavelengths, is 1.6%, 8.8%, 
5.6% and 15.0% for b through e respectively, consistent 
with the photometry error values given in Table [2] The 
largest discrepancies appear in the fake source placed at 
the radial distance of e, as expected, due to its close 
proximity to the star. In particular the J-band peak is 
over-estimated, although the broad-band shape and slope 
of the spectra are retrieved to better than 20% accuracy. 
We note that all molecular features are reproduced, ex- 
cept for the very reddest part of the i^-band, which is 
noisy in any case (cf. Fig. |4|. The critical point with 
respect to interpreting the real spectra of the compan- 
ions is that every minor feature (including, for example, 
the small dent in the J-band peak at about 1250 nm) 
is reproduced with fidelity in the simulated source spec- 
tra. Generally Y and J-band data are worse, because of 
poorer AO performance at those wavelengths. 

7.4. Low- Order Bias Tests 



Because of the possibihty of apparent systematic, but 
broad-band, over, or under, estimates of the peaks es- 
pecially for the radial location of e, as described in the 
last section, we investigated the behavior of injected fake 
sources with perfectly white spectra positioned at vari- 
ous locations in the field of view. This test was designed 
to reveal any field-dependent biases in the spectral ex- 
traction by S4 and to ascertain whether a further minor 
correction to the spectra extracted is necessary. No such 
correction was deemed necessary, as we explain below. 

We injected fake white-spectrum companions at a 
range of intensities (expressed as a fraction of the lo- 
cal speckle intensity) and at a range of radial distances 
from the star (from 15 to 50 pixels, in steps of 5 pixels). 
For each radial location we injected fake sources with a 
range of 0.5% to 5% intensity relative to the local speckle 
brightness, with steps of 0.5%. For each of these 140 fake 
sources, we retrieved the spectrum, as we did with the 
tests in §7.3[ and computed the rms error between the 
injected and retrieved spectra as well as the difference 
in slope between the two. None of the extracted spectra 
had a linear slope that deviates from white by more than 
0.2 rms. To visualize the results of this test, we present 
a surface plot of the rms error versus radial location and 
brightness Fig. [7[ The rms error at the locations of each 
of the four planets is less than 0.4 in all cas es, c onsistent 
with the tests utilizing a T4.5 spectrum (^7.3). No lo- 
cation in the explored parameter space exceeds an rms 
error of 0.75. 

Since no systematic error was detected, and no addi- 
tional correction of spectral shape is required, we con- 
clude that the spectra shown in Fig. [4] are real and as 
accurate as the data permit. Thus, it is appropriate to 
interpret them. 

8. OBSERVATIONAL INTERPRETATION OF THE 
SPECTRA 

Even in comparison with other very low temperature 
objects, such as T-dwarfs, the spectra in Fig. [4] are, first, 
all different from each other and have significant differ- 
ences from known objects. The HR 8799 companions all 
exhibit red colors in the near IR, while the lowest tem- 
perature brown dwarfs maintain relatively blue colors in 
J and iJ-bands. However, the newly discovered exam- 
ples of so-called "Y" dwarfs, though likely much cooler 
than the objects studied in this paper, may b e red as 



shown by [T inney et al.] ( |2012|), [l 
et al.| ( 2012p and especially Kirk; 



Liu et al. (2012LMQ^ley 



Patrick etar] l |2(J12[ ). in 



fact it has been suggested in some of these papers that as 
sensitivity to very low temperature or low surface grav- 
ity objects increases, instruments may need to be built 
to observe at increasingly longer wavelengths, with the J 
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Fig. 6. Spectrum of the T4.5 standard 2MASS J0559-1404 (solid black line, labeled "ground truth;" [Burgasser Kirkpatri^|2006| > 
down-resolved to the Project 1640 resolution and used as a test of the spectral extraction technique. A fake source was injected into tiie 
data with this spectrum at the same radial location and intensity as each of the planets in the HR 8799 system. The retrieved spectrum 
after the S4 extraction is shown as the colored solid line and with the colored dotted lines indicating the 1-cr error bars. Each extraction 
is consistent with the original spectrum in slope and individual features (see This demonstrates the high fidelity of the S4 spectral 
extraction technique. 



TABLE 4 

Statistical Covariance for All Pairs of Spectra 



b and c b and d b and e c and d c and e d and e c and c^ 
0.152 0.015 -0.050 0.045 0.027 0.003 0.849 



^ For all but the first column only points in the if-band are included. The last 
column is the comparison of c between the two epochs of observation (June 
2012 and October 2012). 



and i^-bands becoming so under-luminous that the ob- 
jects of interest fah below detection limits. Indeed, the 
d component exhibits a J—H color of > 2"^, whereas a 
mid-T dwarf has di J—H color of 0.0. However, it is im- 
portant to note that some of the warmer L-dwarfs can 
have colors that are this red in the near-IR. 

8.1. Are the Four Spectra Statistically Different? 

As a very rudimentary initial measurement, we calcu- 
lated the statistical covariance among the six possible 
pairs of four spectra and the two spectra of c, to estab- 
lish how similar or dissimilar they are from each other. 
Except for the comparison of the two epochs of the ob- 
servation of c, no two have a covariance above 0.8 (a 



common measure of a real correlation) and all are below 
0.2, a general benchmark for a significant statistical dif- 
ference. Components b and c have a covariance of 0.152, 
which also indicates significant statistical difference be- 
tween the two objects, though their covariance is higher 
than in the other cases. The calculated values are tab- 
ulated in Table [U For the calculations in which d or e 
were involved, only i^-band data were used, since the Y 
and J spectra are only marginal detections. 

The covariance values suggest that a very weak similar- 
ity exists between b and c (essentially the main hump in 
the i^-band). These values also suggest that the others 
are all significantly different from each other, confirming 
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RMS Spectral Extraction Error 




covariance of these two spectra is —0.41 which is domi- 
nated by the broad slope difference. 
In the bottom- left panel of Fig. [S] w e overlay the spec- 



t 
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Fig. 7. — Surface plot showing the rms error on extracted white 
spectrum sources injected into the data at a range of intensities 
(0.5% to 5%) relative to the local speckle brightness in the unpro- 
cessed data, and a range of radial distances from the star (in pixels, 
ranging from 15 to 50). The positions of the detected planets in 
this parameter space are shown by white crosses and labeled. This 
test demonstrates the high fidelity of the S4 spectral extraction 
technique. 

an inspection by eye of Fig. [4j 

8.2. Comparison with Other Objects 

For comparison purposes, we attempted to find spec- 
tra of other celestial point sources that exhibit at least 
some similarity to those of the component s of HR 8799. 
We pr imarily used the library of spectra in |Rayner et al. 
( 2009 ). Initially we plotted the four spectra overlaid with 
those of the outer planets of the solar system. However, 
the planets do not provide reasonable matches, with the 
sole exception of the e component, as described below. 
We expanded the pool of comparison objects to stars 
and brown dwarfs. In Fig. [S] we show comparisons with 
objects of best match presented for each of the four com- 
panions of HR 8799. In all cases there are discrepancies, 
which clearly shows that these four objects are currently 
unique. Sources b and c seem to exhibit the usual water 
features between F, J and i^-bands, while for d and e we 
only detect the red edge of the water band between J and 
iJ-band. The main methane feature at 1.65 /im seems 
present only in the d and e components, and strongest 
in e. See §8.3| 

In Fig. Is] (top-left), we overlay the spectrum of the 
peculiar L5 brown dwa rf 2MASSJ035523.51+113337.4 



(2M0355; [Faherty et al ||2013J with the b component of 
HR 8799. I'his brown dwart, sometimes referred to with 
a L57 spectral type, seems best explained with low grav- 
ity a nd a possible reddening due to dust (Faherty et al. 
|2013). Although the general shape of the spectrumTof 
zM0355 may be a match, none of the small features are 
reproduced. We believe these to be molecular features 
and they are discussed in §8.3[ 

In the top-right panel of Fig. [8) we compare HR 8799c 
to the T2 dwarf SDSS J125453.90-012247.4. Clearly this 
object is not completely similar to a field T dwarf, as sig- 
nificant discrepancies appear at many wavelength chan- 
nels. The comparison T dwarf is somewhat bluer with 
HR 8799c showing a lower brightness in the F-band, and 
the lack of methane in the i^-band, while the J-band 
peak appears sharper, or "pointier." We note that the 



trum of 2M0355 ( [Faherty et al.|2013D with the d compo- 
nent of HR 8799, as we did in tne top-left panel for the b 
component. Important discrepancies exist in several fea- 
tures around 1500 nm and 166 nm in the d spectrum. 
These are discussed further in ^8^ 

The bottom-right panel of Fig. 8] shows perhaps the 
best match between another object and one of the 
HR 8799 planets, the e component. Overlaid is the spec- 
trum of Saturn, but reddened by a function of and 
normalized to match the overall fiux level. This red- 
dening function is really an attempt to approximate the 
spectrum of Saturn with the Sun's Rayleigh- Jeans black 
body spectral tail removed, as though it were a night-side 
spectrum of Saturn. We are unaware of any night-side 
spectra of the giant planets in this wavelength range, and 
so employed this simplistic approach. We note that there 
is a strong indication of methan e in this object. 

Recently, [ Marley et ah] ( |2012| ) and [Morley et al] ( |2012 ) 
have produced new models that incorporate far more 
complex and realistic cloud and condensate behavior. 
These have been applied to the previously published pho- 
tometry of the HR 8799 system and indicate that clouds 
at varying levels in the atmosphere and with incomplete 
covering fraction of the atmosphere (as is seen in the 
giant planets of our solar system) may explain the red 
behavior of the spectral energy distributions. The new 
data presented here should provide additional constraints 
to these models. 

8.3. Initial Identification of Molecular Features 



Using the database of Freedman et al. (2008) and 
in personal communications with iVl. Marley and D. 
Saumon, we have indicated in Fig. [4] the four main 
sources of opacity that most closely match the features 
we see in these four objects (aside from the water opac- 
ities between the astronomical bandpasses). These in- 
clude ammonia, NH3, in the 1450 to 1550 nm range, 
acetylene, C2H2, in the 1500 to 1550 nm range, methane, 
CH4, redward of 1650 nm, and possibly CO2 in the 
1560 to 1630 nm range, a feature which is commonly 
used in terrestrial atmospheric science (e.g. O'Brie nfc] 
Rayner| [200 2 ) . Given the temperatures and masses of 
these objects as discussed in ^ these seemed the sim- 
plest molecules to consider. In Fig.[4]we indicate the CO2 
feature with a question mark because it is unclear that 
it should be visible at the high temperatures inferred for 
these objects. 

Based on the spectra, and assuming we have identified 
the sources of opacity correctly, we propose the following: 

• b: contains ammonia and/or acetylene as well as 
CO2 but little methane. 

• c: contains ammonia, perhaps some acetylene but 
neither CO2 nor substantial methane. 

• d: contains acetylene, methane and CO2 but am- 
monia is not definitively detected. 

• e: contains methane and acetylene but no ammonia 
or CO2. 
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Spectra of HR 8799 b (top left), c (top right), d (bottom left) and e (bottom right) with comparison objects. 



We emphasize here that these are tentative identifica- 
tions, although the lack of methan e in the b compo nent 
is consistent with the work ini Skemer et al. (2012), Bar- 
an et al. (| 2Qll| a n d Bowler et al. (2010 ). In contrast. 
Barman et al. (201]J) identify CO, not CO2 as we indicate 
Here. Detailed modeling by theorists will be necessary to 
definitively test these identifications. 

Acetylene, has never been convincingly identified in a 
substellar object outside the solar system, and unfortu- 
nately the fundamental band is in the thermal infrared, 
for which we have no spectral coverage. Therefore it is 
difficult to confirm with other observations at this point. 

8.4. Variability and Source of Some Molecules 

As mentioned in ^ there is a weak 2-a difference in 
the c spectrum between June 2012 and October 2012 
that is suggestive of the presence of a weak CO2 feature 
at ~ 1610 nm in October 2012 that is not apparent in 
June 2012. Given the new models with patchy clouds, 
which imply variability in these objects, it is not a leap of 
faith to consider variability of certain molecular features 
in these objects. It is probably unlikely that broad-band 
variability in the primary star (as discussed in Q is the 
cause of such a feature, mainly because the emergent ffux 
from each of these planets is on the order of a few times 
10 ergs cm~^ s~^ (assuming the ~ 800K temperatures 



others have derived) while the incident starlight ranges 
from 3 to 4 orders of magnitude smaller. 

However, the UV ffux from the primary star is consid- 
erable because it is an A5V star — with more than 1000 
times the UV ffux of the Sun. If the Lyman-a line is 
also prominent in HR 8799, there may be larger UV ra- 
diation incident upon even the most dist ant b compo- 
nent c ompared to that upon Jupiter. 



aD 



As [Zahnle et 

(|2009|) demonstrate, UV ffux combined with strong mix- 



ing m the atmosphere can greatly complicate the chem- 
istry represented in the emergent spectra of planets (such 
as Jupiter, in their paper). Acetylene and even hydro- 
gen cyanide (HON) can become abundant. We mention 
HON because it has an absorption feature at 1580 nm, 
comparable to the feature we have tentatively assigned 
to CO2. A known FeH feature exists there as well, but 
has broader extent. 

We will continue to monitor the spectra of these plan- 
ets into the future to assess the reality of any vari- 
ability. We note that for years photometric variability 
of L-dwarfs has be en k nown, and detected in some T- 
dwarfs as well (e.g. [Khandrika et al.|2Q13j ^ and references 
therein). From this point of view, it would not be sur- 
prising to detect variability in any or all of these planets. 

As discussed in ^ and shown in Table [2| c, d, and 
e differ by less than a magnitude in i^-band, but their 
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spectra are not the same. Furthermore, b has some sim- 
ilarity in spectral shape when compared with c, which is 
nearly a magnitude brighter. This includes the peak near 
1350 nm and the general shape of the large feature in H- 
band. In addition, the e component has similarity with a 
reddened spectrum of Saturn. All these facts point to a 
diversity in salient properties of objects for which mass, 
age and metallicity are not the only factors determining 
observables. 

Our results also highlight the power of discerning spec- 
tral differences between planets in multi-planet systems, 
which can provide constraints on the formation loc ations 
for each object. Specifically, Oberg et al. (2011) point 



out that the measured carbon-to-oxygen ratios m exo- 
planet atmospheres may constrain the formation loca- 
tions of such objects, since various snow lines of carbon 
and oxygen-rich ices form at various radial locations from 
the star. 

The results presented here consist of the first compar- 
ative spectroscopic study of multiple planets around a 
star other than our Sun. 

Finally we note that it would be of tremendous value to 
the exoplanet community for any space missions passing 
behind any of the planets of the solar system to take 
night-side thermal spectra of the planets in the near IR, 
to enable direct comparisons with exoplanets. 

9. DISCUSSION: PLANETS ARE DIVERSE 

We have attempted to provide significant new material 
in the understanding of the exosolar system of HR 8799, 
through spectroscopy of all of the known planets. What 
is most striking is that previous studies of the system in- 
dicate that components c, d and e are roughly the same 
near-IR luminosity (within roughly a magnitude). How- 
ever, as we have shown, their spectra are substantially 
different from each other. What causes this remains un- 
known. Whether this is due to small differences in for- 
mation, metallicity differences as a function of orbital 
radius, or evolutionary differences all remain questions. 
On the other hand, it is important to note that these are 
the first spectroscopic observations of multiple planets 
in a planetary system other than our own. Thus it is, 
perhaps, not surprising that four planets that differ by 
only a magnitude or two in brightness exhibit such diver- 
sity in their spectra. Furthermore, similar rapid spectral 
changes happen at the L/T and M/L transitions. 

The spectra of HR 8799 b, c, d and e indicate much 
redder colors than objects with similar spectral features 
(such as methane) currently known. Some authors have 
provided explanations for the extremely red nature of the 
spectra through far higher cloud content than previously 
thought for objects in this mass and temperature range. 
In addition, a cloud covering fractio n smaller than 100% 
may be exhibited in these objects dSkemer et al.] 2012 



., ( I Skemer et al.||2012 

Currie et al.|2011[|Marley et al.|2012[|lVlorley et al.|2012p 

Given the portrait ot the star and other aspects ot the 
system described in ^T] we note, in particular, the ex- 
treme variability of the star in the context of our own 
solar system. The primary star can vary by as m uch as 
8% in t/- b and luminosity over a period of two days ( Zerbi 



et al. 1999). Compared to the Sun, this is a violent en- 



vironment for the stability of any planet's atmosphere. 
The solar constant exhibits v ariability on the order of 
several parts per million (e.g. Froehlich 1987 Fedorov 



2012), and the variability due to the Earth's non-circular 
orbit is on the scale of 6%, but over 6 months, not two 
days ( |Fedorov|[2QT2| . In the case of HR 8799 variabil- 
ity of the star is 3 orders of magnitude greater over a 
hundred times shorter timescale. Furthermore, given the 
A5V spectral class, the planets in this system must be 
irradiated with far larger fluxes in the aerosol-generating 
UV wavelength range. If complex UV photochemsitry 
is present, the planets may exhibit variability in spectral 
features. While our initial detection of such variability in 
CO2 or HCN in the c component is of such low confidence 
that it cannot be claimed definitively, future monitoring 
of the spectra of these objects is warranted. 

10. END NOTE: PROJECT CONTEXT AND SURVEY 

This communication represents the culmination of over 
a decade of work by the authors, including efforts in 
science-grounded instrument conception and design; op- 
tical, mechanical and electrical engineering; development 
of novel techniques for the manipulation and control 
of light from distant stars at the level of A/1000; sys- 
tems engineering and integration; control and data re- 
duction software, which, in this case, comprises several 
large efforts that should be considered "instruments" in 
their own right; software for the identification and spec- 
trum extraction of possible companions — an effort that 
includes expertise from the field of computer vision; ad- 
vanced detector control; and all of the tools of modern 
astronomy brought to bear on the fundamentally difficult 
problem of high-contrast imaging: astrometry, coronag- 
raphy, spectroscopy, photometry and various aspects of 
point source analysis and signal processing. 

This paper marks the beginning of a three-year survey 
to search for and characterize new exoplanetary systems 
around the nearest A and F stars. With the advent of 
this new project and several others beginning operations 
in the next two years, comparative exoplanetary science 
is beyond the initial technical hurdles and can now move 
into a stage of exploring the range of planets extant. The 
techniques will, of course, continue to improve, and an 
exciting research field lies ahead. 
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APPENDIX 

A. THE KLIP SPECKLE SUPPRESSION ALGORITHM 
A.l. Detection 

The extracted data-cubes are first cleaned for remaining spurious bad pixels using median filtering and high-pass 
filtered to remove the contribution of the resid ual atmospheric ha lo. Their radial sc aling and relative registration are 
then calculated using the procedure discussed in Crepp et al. (2011 ) and Pueyo et al. (2012a). For each wavelength from 
Ao = 995 nm, all of the subsequent slices in each cube of the observing sequence are then compressed or stretched and 
registered. In these pre-processed cubes all speckle patterns appear at the same spatial scale and the PSF at Aq is the 
true image of the sky. This yields a series of Nx pre-processed cubes of dimen sions Nx x A/'Exposnres x A^Spaxeis x ^Spaxeis- 



We then proceed to speckle removal using the KLIP algorithm described in Soummer et al. (2012). The target images 
consist of the A/'sxposures PSFs at Aq in the pre-processed cube uxq - We then proceed to pa rtition the im ages into search 
zones, tS, that are equivalent to the optimization zones in LOCI as presented by Lafreniere et al. ( 2007[ ) (radial location 
r, radial width Ar, azimuthal location ^, azimuthal width A^). Note that in contrast to LUCi, lor KLIP these are 
not subtraction zones. For each target image the reference ensemble is chosen as the subset of rescaled slices in the 
same exposure whose wavelength A is such that the image of a putative planet in the compressed/str e tched images 
is sufficiently far from its location at Aq (parameter Ns in Lafreniere et al. ( 2007[ ); [ Crepp et"aL (2011); Pueyo et al. 
( [2012ap ). Finally the A'Exposures speckle-reduced images at Aq are co- added and we proceed to the next wavelength. 

We conduct a Monte-Carlo search over the ensemble of parameters (A^^, Ar, A^, i^KLip) with Ns = 
0.5,0.75,1,1.25,1.5 W (W is the un-occulted PSFs FWHM), Ar = 10,20,30 pixels, ^ = ^ i^' i i^' i^' ^' ^' 
and i^KLiP = 1.25. This approach is reminiscent of the parameter search in Soummer et al. (2011). However, when 
compared to LOCI, it is greatly facilitated by the absence of subtraction in KLIP and by our moderate number of 
reference PSFs. HR 8799 b and c are visible in iY-band reduced images (SNR > 8) and detected using matched 
filtering with un-occulted PSFs, at SNR ~ 3, in J-band reduced images for ^ < ^ and i^xLiP > 10. Since HR 8799 
d and e lie at closer angular separations they are detected in i^-band by matched filtering, at SNR ~ 3, for ^ > | 
(values of 7V^, see [Lafreniere et aL (2007), similar to the ones at the location of HR 8799 b, c) and i^KLiP > 10- These 
findings are sumniarized in the integrated i^-band image shown in Fig. [2] 

A. 2. Spectrum Extraction 

Large biases of the spectral information can arise wh en estimating the sp ectro-photometry of faint companions 
unravelled by an aggressive PSF subtraction routine. In Pueyo et al. (2012a) we identified two sources of potential 
biases in IFS data: self-subtraction from over fitting the companion's signal, and wavelength-to-waveleng th cross 
talk. When using KLIP the former can be calibrated using the forward modeling methodology discussed in Soummer] 



et al. ([2012|, provided that the self-subtraction is not too severe (e.g. provided that it does not radically change 
the niorphology of the PSF). However since our detection pipeline relies on search zones of large radial extent, some 
companion flux is present in the basis-set of principal components. Choosing a small Kklip can, in principle, mitigate 
the wavelength-to- wavelength cross talk, because the contribution of the companion's flux to the reference PSF is most 
likely in the small eigenvalue components. However, this does not fully alleviate this phenomenon. The full details of 
this characterization pipeline will be reported in a subsequent paper focusing on the astrometric characterization of 
the HR 8799 using Project 1640. Below we only outline the main steps of our method. 

Once the rough location (rc,0c) of the companions is known from the detection pipeline, our spectral extraction 
pipeline is composed of the following steps: 

1. Minimization of wavelength-to- wavelength cross talk: We partition the images in such a way that ensures that 
there is no companion signal in the search zones. For the wavelengths Aq < ^Mid^ where Xuid = 1380 nm is 
the central wavelength of the Project 1640 bandpass, we choose search zones azimuthally centered at (j)c over 
the radial interval [rc — NsW^rc — NsW -\- Ar] and only use reference images with A > Aq. For the wavelengths 
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Ao > ^Mid we choose search zones azimuthally centered at (j)c over the radial interval [tc -\- NsW — Ar, Tc -\- NsW] 
and only use reference images with A < Aq. This ensures that there is no companion signal in the search zones but 
dramatically reduces the number of PSFs in the reference ensemble when only using slices from same exposure. 
We alleviate this issue by aggregating slices at all the relevant wavelengths from the full observing sequence in 
a larger reference ensemble. We then proceed through KLIP reduction for these search zones over the range of 
parameters for which the companion has been previously detected. 

Calibration of the self-subtracti on usin^ forward modeling: Using un-occulted PSFs we apply the forward mod- 
eling methodology described in Soummer et al. (2012). We use an un-occulted PSF at each wavelength as a 
template. For each reduction parameter this process produces single wavelength detection maps of width 9 pixels 
around (rc, ^c)- We then derive the spectro-photometry in each channel as value at the peak of an area of the de- 
tection map contained in a small circle of diameter 3 pixels centered at (rc, ^c)- The random spectro-photometric 
uncertainty is estimated as the scatter within the 63% confidence interval of the detection map around this peak. 
The rather large 9x9 spatial extent of the detection maps is chosen as a sanity check to rule out the cases of 
residual speckles contaminating the companion's signal. This is necessary to rule out the pathological case for 
which the tail of a speckle located ~ 2 pixels away from the companion is considered to be companion flux by 
the matched filter. 



3. Estimation of the systematic uncertainty due to the reduction: The previous steps yield a series of spectra 
over a large collection of reduction parameters (A/'^, Ar, A^, Kklip)- Direct inspection of the detection maps 
determine whether or not the companions have been detected for a given set of parameters. For a given sub-set 
of {Ns, Ar, A^) where detection occurs there are three regimes of PC A truncation: i^KLiP < ^KLiP,Min for which 
residual speckles in the neighborhood of the companion are still visible; i^TKLiP > ^KLiP,Max for which detection 
occurs but the morphology of the PSF has been significantly altered by the self-subtraction and thus cannot 
be calibrated by forward modeling; and i^KLiP,Min < ^klip < ^KLiP,Max which are of scientific interest. We 
analyze our extracted spectra in conjunction with single wavelength detection maps in order to establish the 
subset of parameters over which each planet is detected. The spectra are then derived as the mean of this subset 
of spectra and their uncertainty due to the reduction is derived as the scatter of this ensemble — for HR 8799 b 
and c this subset is rather large (~ 100 spectra) while is is smaller for HR 8799 d and e (~ 10 spectra). 

3. Spectral calibration: The resultant spectra from the previous steps are expressed in fractional units of un- 
occulted core intensity of the primary star (e.g units o f con trast). We derive the spectrum of HR 8799 using 
the A5V teniplate standard from the Pickles library ( |Pickles||l99 8) normalized at mn = 5.28, or Mh = 2.30 



(Cutri et al. 2003| ) and binned at the Project 1640 resolution' The final spectra are obtained by multiplying 
the results from the previous step with this stellar spectrum. Note t hat this procedure for spectral cal ibration 
is slig htly different than in earlier Project 1640 results presented in Hinkley et al. (2010) and [Roberts et al. 
(2012). Those relied on the ratio of aperture photometry estimate s of bo th t he FSl^ cores and com panion. As a 



sanity check we conducted an analysis analogous to Hinkley et al.| (2010) and Roberts et al.] (2012): we used the 
relevant reduced images selected in step 2 and conducted step 3. By projecting un-occulted r]st s on the relevant 
principal components and using aperture photometry to calibrate the self-subtraction, we computed the spectral 
calibration using aperture photometry of un-occulted images. The two yielded results that are consistent within 
our uncertainties. 



B. THE S4 SPECKLE SUPPRESSION ALGORITHM 
B.l. Detection 



The S4 algorithm for post-processing speckle suppression (^Fe rgus et al.|2013" ), based on principal component analysis 

it he 



(PC A), has not yet been published. Thus we provide some details of it here. 

S4 takes as input a 4D data block of dimension A/pixeis x ^pixeisx ^ ^ ^Exposures which has been pre-processed 
as follows: (i) application of a 2 x 2 median filter to each band/exposure to remove dead pixels; (ii) spatial alignment 
of all bands and exposures relative to one another and (iii) spatial centering so that the star lies precisely (within 0.1 
pixel) at (7Vpixeisy/2,iVpixeisx/2). 

In S4 the data at each spatial position is decomposed into a speckle component and a companion component, whose 
sum reconstructs the original data to the limit of Gaussian noise. To model the speckle component effectiv ely, which 
evolves radially with wavelength, this decomposition is performed in a polar reference frame, as shown in Fig. |B1| (left). 
To examine a location at radius d and angle from the center, an annular region of width at a radius d from is 
transformed into a region of size Q x R x Nx x ^Exposures , where B = 27rd to ensure the region is well-sampled. This 
region is divided into a test zone around the location 6 (of size S6) and a training zone of all other angles, as shown 



in red and green respectively in Fig. Bl (left). We use the training zone to build a model of the speckles that is then 
applied to the test region, decomposmg it into speckle and companion components. 

We treat the speckles in both the training and test zones as being independent over angle and exposure n, justified 
by limited angular extent of the speckles and their variation due to atmospheric turbulence and instrument flexure 
or temperature variations between exposures. We thus assume t he st ructure of the speckles to be confined to a joint 



radius- wavelength space (of dimension Ai?), as illustrated in Fig. Bl (right). This joint space is modeled using PCA, 
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Fig. Al. — Comparison of spectra extracted with the KLIP algorithm and the S4 algorithm (in red). Spectra are normalized to 1 and 
shifted by a constant for comparison. Points that are weaker than a 2-cr detection are excluded in the normalization and all points in Y 
and J bands are excluded for planets d and e. Only a few of the 94 points plotted for each algorithm are discrepant by more than 2-cr, 
indicating good agreement between the two techniques 
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Fig. B1. — Schematic of the S4 Algorithm. Left: Mean over wavelength A and exposures for an example star, FU Ori, whose 
companion is obvious, along with a polar representation of the annular region at radius d. Green and red regions shown training and test 
zones, respectively. Right top: The PCA basis W, computed from the training zone, showing the diagonal structure of the speckles in 
the joint radius- wavelength space. Right bottom: Reconstructions using the PCA model with K = 30 components, (a): Original data 
slice, containing only speckles, (b): Reconstruction using PCA model with k = 30. (c): Error residual. Note the lack of structure, (d): 
Companion model, which has low correlation with (c). Next row: (e): Slice containing companion, (f): Reconstruction of PCA model, 
(g): The error residual shows clear structure associated with the companion, i.e. the PCA speckle model cannot reconstruct the companion 
signal, (h): Companion model, which has a high correlation with (g). 



which approximates each Ai?-dimensional shce of the data at angle cube n as a hnear combination of K orthonormal 
basis vectors W = [wi^ . . . , wk] {K being a user-defined parameter). These basis vectors capture the majority of the 
variance in the radius-wavelength space and are computed by performing an eigendecomposition of the covariance 
matrix built from the training region, reshaped into a matrix of XR dimensions by (O — ^^)A^Exposiires samples. We 
then use the PCA basis W to infer the speckle component of the test zone. W is visualized in Fig. |B1| (upper right). 

In detection, we first use the PCA basis W to fit the data in the test zone and then correlate the residual error with 
a fixed companion model P (previously obtained by calibration with a point white light source). The response for 
position (^, d) is stored in a correlation map. The intuition is that the PCA basis W will effectively model the speckles, 
but not the companion, thus the residual should only contain the companion signal which will respond strongly to the 
matched filter P (see Fig. Bl, lower right). The location of the test zone is systematically moved across all angles at 
a given radius, with the PCA basis being recomputed at each location, because the training zone also changes. The 
process is then repeated for a new radius (i, until all spatial locations in the visual field have been examined. 
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Finally, a normalization is performed on the correlation map that compensates for the varition in flux with radius 
in the original data. The resulting normalized map is then converted back to the Cartesian coordinate frame and is 
the output of the detection algorithm. (See Fig. [3| and |Bl[ c) for examples.) 

B.2. Spectrum Extraction 

Promising peaks in the normalized detection map are then selected for spectrum extraction. This follows the same 
overall modeling approach as detection, except that the spectrum of the companion model P is no longer fixed to be 
white. For a peak at location (d, ^), let the observed data for exposure n be yn (represented as XR dimensional vector). 
We now must estimate both the spectrum of the planet m (an Nx dimensional vector) and the K dimensional PCA 
coefficients Zn for exposure n that reconstruct i/n- Using a Gaussian noise model, this is equivalent to minimizing the 
following convex objective: 

posures 

J2 WVn -{Wz„ + Pm)\\l (Bl) 

n=l 

where Wzn is the speckle component for exposure n and Pm is the companion component (constant across all expo- 
sures). We impose a non-negativity constraint on m, since negative spectra are not physically plausible. Minimization 
is performed using standard optimization software in Matlab. The spectra shown in Fig. |4] are the resulting m vectors 
for the 4 different planet locations. An estimate of the uncertainty is obtained by measuring the variance in m over 
different settings of the model parameters, K. In the case of the data presented here, the number of principal compo- 
nents was varied from 50 to 500 with steps of 25 and the angular fitting size, O, was 3, 5 and 7 pixels. The error bars 
were calculated using principal components of 150, 175 and 200 and of 5 and 7 pixels, settings which yielded the 
best results. For October, due to significantly worse observing conditions, more principal components were needed for 
detection and spectrum extraction: 250, 300 and 350 with 6 of 5 and 7 pixels. 
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